Project - Concrete Strength Prediction¶
Problem statement¶
Objective¶
To predict the concrete strength using the data available in file "concrete.csv". Apply feature engineering and model tuning to obtain 85% to 95% accuracy.
Resources Available¶
The data for this project is available in file https://archive.ics.uci.edu/ml/machine-learning-databases/concrete/compressive/. The same has been shared along with the course content.
Steps and Tasks:¶
Exploratory Data Quality Report Reflecting the Following:
- Univariate analysis – data types and description of the independent attributes which should include (name, range of values observed, central values (mean and median), standard deviation and quartiles, analysis of the body of distributions/tails, missing values, outliers, duplicates(10 Marks)
- Bi-variate analysis between the predictor variables and between the predictor variables and target column. Comment on your findings in terms of their relationship and degree of relation if any. Visualize the analysis using boxplots and pair plots, histograms, or density curves. (10 marks)
- Feature Engineering techniques(10 marks)
- Identify opportunities (if any) to extract new features from existing features, drop a feature(if required) Hint: Feature Extraction, for example, consider a dataset with two features length and breadth. From this, we can extract a new feature Area which would be length * breadth.
- Get the data model ready and do a train test split.
- Decide on the complexity of the model, should it be a simple linear model in terms of parameters or would a quadratic or higher degree.
Creating the Model and Tuning It:¶
Algorithms that you think will be suitable for this project. Use Kfold Cross-Validation to evaluate model performance. Use appropriate metrics and make a DataFrame to compare models w.r.t their metrics. (at least 3 algorithms, one bagging and one boosting based algorithms have to be there). (15 marks) Techniques employed to squeeze that extra performance out of the model without making it overfit. Use Grid Search or Random Search on any of the two models used above. Make a DataFrame to compare models after hyperparameter tuning and their metrics as above. (15 marks)
Attribute Information:¶
Given are the variable name, variable type, the measurement unit, and a brief description. The concrete compressive strength is the regression problem. The order of this listing corresponds to the order of numerals along the rows of the database.
Name Data Type Measurement Description
- Cement (cement) quantitative kg in a m3 mixture Input Variable
- Blast Furnace Slag (slag) quantitative kg in a m3 mixture Input Variable
- Fly Ash (ash) quantitative kg in a m3 mixture Input Variable
- Water(water) quantitative kg in a m3 mixture Input Variable
- Superplasticizer (superplastic) quantitative kg in a m3 mixture Input Variable
- Coarse Aggregate (coarseagg) quantitative kg in a m3 mixture Input Variable
- Fine Aggregate (fineagg) quantitative kg in a m3 mixture Input Variable
- Age(age) quantitative Day (1~365) Input Variable
- Concrete compressive strength(strength) quantitative MPa Output Variable
1. Import the necessary packages¶
In [1]:
#Load Libraries
import pandas as pd
import numpy as np
import seaborn as sns
import datetime
import matplotlib.pyplot as plt
%matplotlib inline
In [2]:
import warnings
warnings.filterwarnings('ignore')
Read the dataset into dataframes¶
In [3]:
concrete_df = pd.read_csv("concrete.csv")
concrete_df.head(10)
Out[3]:
In [4]:
concrete_df.tail(10)
Out[4]:
Apply info, shape, describe, and find the number of missing values in the data¶
In [5]:
# Check to see the number of records in the dataset
concrete_df.shape
Out[5]:
Observation: The dataset comprises of 45211 rows and 17 columns¶
In [6]:
concrete_df.info()
In [7]:
concrete_df.dtypes
Out[7]:
In [8]:
# Check to see if data has any missing values
concrete_df.isnull().any()
Out[8]:
Observation: No Null Values are present¶
In [9]:
#Analyze the distribution of the dataset
concrete_df.describe(include = 'all')
Out[9]:
In [10]:
concrete_df.describe(include = 'all').transpose()
Out[10]:
In [11]:
print(concrete_df.columns)
In [12]:
num_cols = ['cement', 'slag', 'ash', 'water', 'superplastic', 'coarseagg', 'fineagg', 'age','strength']
In [13]:
from matplotlib.cbook import boxplot_stats
# Central tendency distribution
i = 1
for col in num_cols:
print('Attribute: ',col)
print(boxplot_stats(concrete_df[col]))
plt.figure(figsize=(10,8))
sns.boxplot(concrete_df[col])
plt.show();
i = i+1
print('-------------------------------------------------------------------------')
In [14]:
# Skewness distribution
i = 1
for col in num_cols:
plt.figure(figsize=(10,8))
sns.distplot(concrete_df[col])
i = i+1
In [15]:
# Check the unique values in each column of the dataframe
concrete_df.nunique()
Out[15]:
In [16]:
import pandas_profiling
concrete_df.profile_report()
Out[16]:
In [17]:
from pandas_profiling import ProfileReport
profile = ProfileReport(concrete_df, title="Pandas Profiling Report")
profile.to_widgets()
Distribution of Target variable¶
In [18]:
concrete_df.strength.nunique()
Out[18]:
In [19]:
concrete_df.strength.value_counts()
Out[19]:
In [20]:
# Plot the distribution of the target attribute
plt.figure(figsize=(8,6))
sns.distplot(concrete_df['strength']);
plt.axvline(concrete_df['strength'].mean(), linestyle="dashed", color="g", label='mean', linewidth=2)
plt.legend(loc="best")
plt.show()
Bivariate Analysis¶
In [21]:
for col in num_cols:
plt.figure(figsize=(12,12))
ax = sns.regplot(x=col,y=concrete_df['strength'], data=concrete_df);
Multivariate analysis¶
In [22]:
sns.pairplot(data = concrete_df,diag_kind ='kde');
In [23]:
# Correlation matrix
corr = concrete_df.corr()
corr
Out[23]:
In [24]:
# Heatmap
sns.set(font_scale=1.15)
fig,ax=plt.subplots(figsize=(16,16))
sns.heatmap(corr, cmap='GnBu',annot=True,linewidths=0.01,center=0,linecolor='white',square=True)
plt.title('Correlation between attributes',fontsize=16)
ax.tick_params(labelsize=14)
- We can observe a high positive correlation between strength and cement
- Also, age and superplastic are the other two factors influencing strength.
- A strong negative correlation exists between superplastic and water
- Positive correlations exist betweensSuperplastic and Fly Ash, Fine Aggregate
In [25]:
# Identify duplicates records in the data
dupes = concrete_df.duplicated()
sum(dupes)
Out[25]:
In [26]:
duplicate = concrete_df[concrete_df.duplicated()]
print("Duplicate Rows :")
duplicate
Out[26]:
In [27]:
# Removing Duplicates
concrete_df = concrete_df.drop_duplicates()
In [28]:
dupes = concrete_df.duplicated()
sum(dupes)
Out[28]:
In [29]:
# Check to see the number of records in the dataset
concrete_df.shape
Out[29]:
In [30]:
Q1 = concrete_df.quantile(0.25)
Q3 = concrete_df.quantile(0.75)
IQR = Q3 - Q1
print(IQR)
In [31]:
np.where((concrete_df < (Q1 - 1.5 * IQR)) | (concrete_df > (Q3 + 1.5 * IQR)))
Out[31]:
In [32]:
# Check for outliers in dataset using boxplot
fig,ax = plt.subplots(figsize=(35,20))
sns.boxplot(data=concrete_df,width=0.5,palette='GnBu',fliersize=7);
In [33]:
concrete_df_new = concrete_df.copy()
In [34]:
# Central tendency distribution
i = 1
for col in num_cols:
print('Attribute: ',col)
Q1 = concrete_df[col].quantile(q=0.25)
Q3 = concrete_df[col].quantile(q=0.75)
lower_limit = Q1-1.5*(Q3-Q1)
upper_limit = Q3+1.5*(Q3-Q1)
print('Lower Limit: ', lower_limit)
print('Upper Limit: ', upper_limit)
outliers = len(concrete_df.loc[(concrete_df[col] < (lower_limit)) | (concrete_df[col] > (upper_limit))])
print("Total number of Outliers:",outliers)
i = i+1
print('-------------------------------------------------------------------------')
- The attributes 'cement','ash' and 'coarseagg' has no outliers.
- Attributes 'superplastic' and 'strength' has mild outliers/outliers very close to the upper whisker of their corresponding boxplots and as such those outliers will not be treated and left as is.
- Same is the case with the lower end outliers for the attribute 'water' as they lie quite close to lower whisker of the boxplot and hence will not be treated.
- The attribute 'Age' has a varied range, different observations in age column are meant for allowing concrete to undergo rest for a wide range of observations and hence will not undergo outlier treatment.
- All the remaining outliers in the other attributes will be treated by capping them with min/max replacement.
In [35]:
out_cols = ['slag','water', 'fineagg']
In [36]:
#Replacing the outliers by Min/Max Replacement:
for col in out_cols:
q1 = concrete_df_new[col].quantile(0.25)
q3 = concrete_df_new[col].quantile(0.75)
iqr = q3 - q1
min_val = q1-1.5*iqr
max_val = q3+1.5*iqr
concrete_df_new[col].loc[concrete_df_new[col] > max_val] = max_val
In [37]:
concrete_df_new.shape, concrete_df.shape
Out[37]:
In [38]:
# Check for outliers in dataset using boxplot
fig,ax = plt.subplots(figsize=(35,20))
sns.boxplot(data=concrete_df_new,width=0.5,palette='BuPu',fliersize=7);
In [39]:
sns.pairplot(data = concrete_df_new,diag_kind ='kde');
Determining feature importance using pearson correlation matrix analysis¶
In [40]:
# Correlation matrix
corr = concrete_df_new.corr()
corr
Out[40]:
In [41]:
# Heatmap
sns.set(font_scale=1.15)
fig,ax=plt.subplots(figsize=(16,16))
sns.heatmap(corr, cmap='BuPu',annot=True,linewidths=0.01,center=0,linecolor='white',square=True)
plt.title('Correlation between attributes',fontsize=16)
ax.tick_params(labelsize=14)
In [42]:
target_corr = abs(corr['strength'])
#Selecting highly correlated features
rel_features = target_corr[target_corr>0.30]
rel_features
Out[42]:
Relationship between indepndendent variables¶
In [43]:
print(concrete_df_new[["cement","superplastic"]].corr())
In [44]:
print(concrete_df_new[["age","cement"]].corr())
- From the pairplot and heatmap above it can be seen that the attributes fineagg,coarseagg, ash and slag have a very week relationship with the target attribute strength
- Cement,water,age,superplastic attrbiutes are more reliable in predicting strength.
Model Building with Feature Selection/Elimination¶
In [45]:
#Scaling the dataset
from sklearn.preprocessing import StandardScaler
from sklearn.metrics import mean_squared_error, r2_score
from sklearn.model_selection import KFold, cross_val_score
from sklearn import metrics
import numpy as np
std_scale = StandardScaler()
concrete_df_scaled = pd.DataFrame(std_scale.fit_transform(concrete_df_new),columns = concrete_df_new.columns)
In [46]:
# Dropping the least significant features
X = concrete_df_scaled.drop(['slag','ash','strength'], axis=1)
#X = concrete_df_scaled.drop(['slag','ash','strength','fineagg','coarseagg'], axis=1)
Y = concrete_df_scaled[['strength']]
In [47]:
# Split into training and test set
from sklearn.model_selection import train_test_split
X_train, X_test, y_train, y_test = train_test_split(X, Y, test_size=0.3,random_state=1)
Linear Regression with K-fold cross validation¶
In [48]:
from sklearn.metrics import confusion_matrix, recall_score, precision_score, f1_score, roc_auc_score,accuracy_score
from sklearn.linear_model import LinearRegression
from sklearn.model_selection import cross_val_score,cross_val_predict
from sklearn.metrics import mean_squared_error, r2_score,mean_absolute_error
linreg = LinearRegression()
In [49]:
linreg_scores = cross_val_score(linreg, X_train, y_train, cv=KFold(10, random_state = 1,shuffle = True),scoring='neg_mean_squared_error')
linreg.fit(X_train, y_train)
print(linreg.coef_)
print()
# Calculate the score of Linear Regression
print('Training score :', linreg.score(X_train, y_train))
print('Testing score :', linreg.score(X_test, y_test))
# Calculate average RMSE
print('RMSE:',np.sqrt(-(linreg_scores)).mean())
In [50]:
#Store the results for each model in a dataframe for final comparison
results_df = pd.DataFrame({'Method':['Logistic Regression'],'Training Score': linreg.score(X_train, y_train).round(3),'Testing score': linreg.score(X_test, y_test).round(3),
'RMSE': np.sqrt(-(linreg_scores)).mean().round(3)})
results_df
Out[50]:
Ridge Regression with K-fold cross validation¶
In [51]:
from sklearn.linear_model import Ridge
ridge = Ridge(alpha = 0.01,random_state=1)
ridge_scores = cross_val_score(ridge, X_train, y_train, cv=KFold(10, random_state = 1,shuffle = True), scoring='neg_mean_squared_error')
ridge.fit(X_train,y_train)
print ("Ridge model:", (ridge.coef_))
print()
# Calculate the score of Ridge Regressor
print('Training Score:',ridge.score(X_train, y_train))
print('Testing Score:',ridge.score(X_test, y_test))
# calculate the average RMSE
print('RMSE:',np.sqrt(-(ridge_scores)).mean())
In [52]:
#Store the accuracy results for each model in a dataframe for final comparison
r_df = pd.DataFrame({'Method':['Ridge Regression'],'Training Score': ridge.score(X_train, y_train).round(3),'Testing score': ridge.score(X_test, y_test).round(3),
'RMSE': np.sqrt(-(ridge_scores)).mean().round(3)})
results_df = results_df.append(r_df)
results_df
Out[52]:
Lasso Regression with K-fold cross validation¶
In [53]:
from sklearn.linear_model import Lasso
lasso = Lasso(alpha = 0.01,random_state=1)
lasso_scores = cross_val_score(lasso, X_train, y_train, cv=KFold(10, random_state = 1,shuffle = True), scoring='neg_mean_squared_error')
lasso.fit(X_train,y_train)
print ("Lasso model:", (lasso.coef_))
print()
# Calculate the score of Lasso Regressor
print('Training Score:',lasso.score(X_train, y_train))
print('Testing Score:',lasso.score(X_test, y_test))
# calculate the average RMSE
print('RMSE:',np.sqrt(-(lasso_scores)).mean())
In [54]:
#Store the accuracy results for each model in a dataframe for final comparison
r_df = pd.DataFrame({'Method':['Lasso Regression'],'Training Score': lasso.score(X_train, y_train).round(3),'Testing score': lasso.score(X_test, y_test).round(3),
'RMSE': np.sqrt(-(lasso_scores)).mean().round(3)})
results_df = results_df.append(r_df)
results_df
Out[54]:
SVM Regression with K-fold cross validation¶
In [55]:
from sklearn.svm import SVR
svr = SVR()
svr_scores = cross_val_score(svr, X_train, y_train, cv=KFold(10, random_state = 1,shuffle = True), scoring='neg_mean_squared_error')
svr.fit(X_train, y_train)
Out[55]:
In [56]:
# Calculate the score of SVM Regressor
print('Training Score:',svr.score(X_train, y_train))
print('Testing Score:',svr.score(X_test, y_test))
# calculate the average RMSE
print('RMSE:',np.sqrt(-(svr_scores)).mean())
In [57]:
#Store the accuracy results for each model in a dataframe for final comparison
r_df = pd.DataFrame({'Method':['SVM Regression'],'Training Score': svr.score(X_train, y_train).round(3),'Testing score': svr.score(X_test, y_test).round(3),
'RMSE': np.sqrt(-(svr_scores)).mean().round(3)})
results_df = results_df.append(r_df)
results_df
Out[57]:
Decision Tree with K-fold cross validation¶
In [58]:
from sklearn.tree import DecisionTreeRegressor
dtree = DecisionTreeRegressor(random_state=1)
dtree_scores = cross_val_score(dtree, X_train, y_train, cv=KFold(10, random_state = 1,shuffle = True), scoring='neg_mean_squared_error')
dtree.fit(X_train, y_train)
Out[58]:
In [59]:
#Making the prediction
y_predict = dtree.predict(X_test)
In [60]:
print (pd.DataFrame(dtree.feature_importances_, columns = ["Imp"], index = X_train.columns))
In [61]:
# Calculate the score of Decision Tree Regressor
print('Training Score:',dtree.score(X_train,y_train))
print('Testing Score:',dtree.score(X_test,y_test))
# calculate the average RMSE
print('RMSE:',np.sqrt(-(dtree_scores)).mean())
In [62]:
#Store the accuracy results for each model in a dataframe for final comparison
r_df = pd.DataFrame({'Method':['Decision Tree'],'Training Score': dtree.score(X_train, y_train).round(3),'Testing score': dtree.score(X_test, y_test).round(3),
'RMSE': np.sqrt(-(dtree_scores)).mean().round(3)})
results_df = results_df.append(r_df)
results_df
Out[62]:
Though the testing score has improved, there is quite a drop in accuracy from training score to testing score indicating this model is highly overfit and needs to be pruned.
Pruned Decision Tree with K-fold cross validation¶
In [63]:
dtree_pruned = DecisionTreeRegressor(max_depth=6,random_state=1)
dtree_pruned_scores = cross_val_score(dtree_pruned, X_train, y_train, cv=KFold(10, random_state = 1,shuffle = True), scoring='neg_mean_squared_error')
dtree_pruned.fit(X_train, y_train)
Out[63]:
In [64]:
#Making the prediction
y_predict = dtree_pruned.predict(X_test)
In [65]:
# Calculate the score of Decision Tree Regressor
print('Training Score:',dtree_pruned.score(X_train,y_train))
print('Testing Score:',dtree_pruned.score(X_test,y_test))
# calculate the average RMSE
print('RMSE:',np.sqrt(-(dtree_pruned_scores)).mean())
In [66]:
#Store the accuracy results for each model in a dataframe for final comparison
r_df = pd.DataFrame({'Method':['Decision Tree with Pruning'],'Training Score': dtree_pruned.score(X_train, y_train).round(3),'Testing score': dtree_pruned.score(X_test, y_test).round(3),
'RMSE': np.sqrt(-(dtree_pruned_scores)).mean().round(3)})
results_df = results_df.append(r_df)
results_df
Out[66]:
Random Forest with K-fold cross validation¶
In [67]:
from sklearn.ensemble import RandomForestRegressor
rforest = RandomForestRegressor(random_state=1)
rforest_scores = cross_val_score(rforest, X_train, y_train, cv=KFold(10, random_state = 1,shuffle = True), scoring='neg_mean_squared_error')
rforest = rforest.fit(X_train, y_train)
In [68]:
y_predict = rforest.predict(X_test)
In [69]:
# Calculate the score of Random Forest Regressor
print('Training Score:',rforest.score(X_train,y_train))
print('Testing Score:',rforest.score(X_test,y_test))
# calculate the average RMSE
print('RMSE:',np.sqrt(-(rforest_scores)).mean())
Similar is the case with Random Forest regressor. Though the testing score has improved, there is quite a drop in accuracy from training score to testing score indicating this model is highly overfit and needs to be pruned.
In [70]:
#Store the accuracy results for each model in a dataframe for final comparison
r_df = pd.DataFrame({'Method':['Random Forest'],'Training Score': rforest.score(X_train, y_train).round(3),'Testing score': rforest.score(X_test, y_test).round(3),
'RMSE': np.sqrt(-(rforest_scores)).mean().round(3)})
results_df = results_df.append(r_df)
results_df
Out[70]:
Pruned Random Forest with K-fold cross validation¶
In [71]:
reg_rforest = RandomForestRegressor(n_estimators = 50,max_depth = 7,min_samples_leaf=5,random_state=1)
reg_rforest_scores = cross_val_score(reg_rforest, X_train, y_train, cv=KFold(10, random_state = 1,shuffle = True), scoring='neg_mean_squared_error')
reg_rforest = reg_rforest.fit(X_train, y_train)
In [72]:
y_predict = reg_rforest.predict(X_test)
In [73]:
# Calculate the score of Random Forest Regressor
print('Training Score:',reg_rforest.score(X_train,y_train))
print('Testing Score:',reg_rforest.score(X_test,y_test))
# calculate the average RMSE
print('RMSE:',np.sqrt(-(reg_rforest_scores)).mean())
In [74]:
#Store the accuracy results for each model in a dataframe for final comparison
r_df = pd.DataFrame({'Method':['Random Forest with Pruning'],'Training Score': reg_rforest.score(X_train, y_train).round(3),'Testing score': reg_rforest.score(X_test, y_test).round(3),
'RMSE': np.sqrt(-(reg_rforest_scores)).mean().round(3)})
results_df = results_df.append(r_df)
results_df
Out[74]:
Bagging with K-fold cross validation¶
In [75]:
from sklearn.ensemble import BaggingRegressor
bag = BaggingRegressor(random_state=1)
bag_scores = cross_val_score(bag, X_train, y_train, cv=KFold(10, random_state = 1,shuffle = True), scoring='neg_mean_squared_error')
bag = bag.fit(X_train, y_train)
In [76]:
y_predict = bag.predict(X_test)
In [77]:
# Calculate the score of Bagging Regressor
print('Training Score:',bag.score(X_train,y_train))
print('Testing Score:',bag.score(X_test,y_test))
# calculate the average RMSE
print('RMSE:',np.sqrt(-(bag_scores)).mean())
In [78]:
#Store the accuracy results for each model in a dataframe for final comparison
r_df = pd.DataFrame({'Method':['Bagging'],'Training Score': bag.score(X_train, y_train).round(3),'Testing score': bag.score(X_test, y_test).round(3),
'RMSE': np.sqrt(-(bag_scores)).mean().round(3)})
results_df = results_df.append(r_df)
results_df
Out[78]:
Bagging with SVR as base estimator¶
In [79]:
bag_svr = BaggingRegressor(base_estimator=SVR(),random_state=1)
bag_svr_scores = cross_val_score(bag_svr, X_train, y_train, cv=KFold(10, random_state = 1,shuffle = True), scoring='neg_mean_squared_error')
bag_svr.fit(X_train, y_train)
Out[79]:
In [80]:
y_predict = bag.predict(X_test)
In [81]:
# Calculate the score of Bagging Regressor
print('Training Score:',bag_svr.score(X_train,y_train))
print('Testing Score:',bag_svr.score(X_test,y_test))
# calculate the average RMSE
print('RMSE:',np.sqrt(-(bag_svr_scores)).mean())
In [82]:
#Store the accuracy results for each model in a dataframe for final comparison
r_df = pd.DataFrame({'Method':['Bagging with SVR'],'Training Score': bag_svr.score(X_train, y_train).round(3),'Testing score': bag_svr.score(X_test, y_test).round(3),
'RMSE': np.sqrt(-(bag_svr_scores)).mean().round(3)})
results_df = results_df.append(r_df)
results_df
Out[82]:
Gradient Boost with K-fold cross validation¶
In [83]:
from sklearn.ensemble import GradientBoostingRegressor
gboost = GradientBoostingRegressor(random_state=1)
gboost_scores = cross_val_score(gboost, X_train, y_train, cv=KFold(10, random_state = 1,shuffle = True), scoring='neg_mean_squared_error')
gboost = gboost.fit(X_train, y_train)
In [84]:
y_predict = gboost.predict(X_test)
In [85]:
# Calculate the score of Gradient Boost Regressor
print('Training Score:',gboost.score(X_train,y_train))
print('Testing Score:',gboost.score(X_test,y_test))
# calculate the average RMSE
print('RMSE:',np.sqrt(-(gboost_scores)).mean())
In [86]:
#Store the accuracy results for each model in a dataframe for final comparison
r_df = pd.DataFrame({'Method':['Gradient Boost'],'Training Score': gboost.score(X_train, y_train).round(3),'Testing score': gboost.score(X_test, y_test).round(3),
'RMSE': np.sqrt(-(gboost_scores)).mean().round(3)})
results_df = results_df.append(r_df)
results_df
Out[86]:
XGBoost with K-fold cross validation¶
In [87]:
import xgboost as xgb
from xgboost import XGBRegressor
xgboost = XGBRegressor(random_state=1)
xgboost_scores = cross_val_score(xgboost, X_train, y_train, cv=KFold(10, random_state = 1,shuffle = True), scoring='neg_mean_squared_error')
xgboost = xgboost.fit(X_train, y_train)
In [88]:
y_predict = xgboost.predict(X_test)
In [89]:
# Calculate the score of XGBoost Regressor
print('Training Score:',xgboost.score(X_train,y_train))
print('Testing Score:',xgboost.score(X_test,y_test))
# calculate the average RMSE
print('RMSE:',np.sqrt(-(xgboost_scores)).mean())
In [90]:
#Store the accuracy results for each model in a dataframe for final comparison
r_df = pd.DataFrame({'Method':['XGBoost'],'Training Score': xgboost.score(X_train, y_train).round(3),'Testing score': xgboost.score(X_test, y_test).round(3),
'RMSE': np.sqrt(-(xgboost_scores)).mean().round(3)})
results_df = results_df.append(r_df)
results_df
Out[90]:
Random Forest, Gradient Boost and XGBoost Regressors seem to be the best suited models with higher accuracy and lower RMSE scores. Let's tune the models and check for model performance at 95% confidence level.
Model Tuning¶
Hyperparameter Tuning - GridSearchCV - Random Forest¶
In [91]:
from pprint import pprint
print('Parameters currently in use:\n')
pprint(rforest.get_params())
In [92]:
# Create the GridSearch estimator along with a parameter object containing the values to adjust
from sklearn.model_selection import GridSearchCV
n_estimators = [200, 500,1000,1800,2000]
max_features = ['auto','sqrt']
max_depth = [8,10,20,30,40,50]
min_samples_split = [2, 5, 10]
min_samples_leaf = [1, 2, 5, 10]
param_grid = dict(n_estimators = n_estimators, max_depth = max_depth,
min_samples_split = min_samples_split, max_features = max_features,
min_samples_leaf = min_samples_leaf)
In [93]:
# Create grid search using 5-fold cross validation
clf = GridSearchCV(rforest, param_grid, cv=5,n_jobs = -1,verbose=2)
In [94]:
# Fit grid search
best_model = clf.fit(X_train, y_train)
In [95]:
# View best hyperparameters
training_score = best_model.best_score_
print(best_model.best_score_)
print(best_model.best_params_)
# Actual model object fit with those best parameters
# Shows default parameters that we did not specify
print(best_model.best_estimator_)
In [96]:
# Apply best params to fit this model and see training scores
#clfg = RandomForestRegressor(random_state=1, max_depth=10, n_estimators=2000, max_features='sqrt',min_samples_leaf=1, min_samples_split=2)
clfg = best_model.best_estimator_
print("R Squared:",cross_val_score(clfg, X_train, y_train.values.ravel(), cv=10, scoring='r2').mean())
In [97]:
clfg = clfg.fit(X_train, y_train.values.ravel())
In [98]:
test_RMSE = np.sqrt(mean_squared_error(y_test, clfg.predict(X_test)))
print("Test RMSE: ", np.sqrt(mean_squared_error(y_test, clfg.predict(X_test))))
test_score = r2_score(y_test, clfg.predict(X_test))
print("Test R^2: ", r2_score(y_test, clfg.predict(X_test)))
In [99]:
#Store the accuracy results for each model in a dataframe for final comparison
r_df = pd.DataFrame({'Method':['Random Forest with GridSearch CV'],'Training Score': training_score.round(3),'Testing score': test_score.round(3),
'RMSE': test_RMSE})
results_df = results_df.append(r_df)
results_df
Out[99]:
Hyperparameter Tuning - RandomizedSearchCV - Random Forest¶
In [100]:
from sklearn.model_selection import RandomizedSearchCV
n_estimators = [200, 500,1000,1800,2000]
max_features = ['auto','sqrt']
max_depth = [8,10,20,30,40,50]
min_samples_split = [2, 5, 10]
min_samples_leaf = [1, 2, 5, 10]
param_grid = dict(n_estimators = n_estimators, max_depth = max_depth,
min_samples_split = min_samples_split, max_features = max_features,
min_samples_leaf = min_samples_leaf)
In [101]:
# Create grid search using 5-fold cross validation
clf = RandomizedSearchCV(rforest, param_grid, cv=5, n_iter = 100,n_jobs = -1,verbose=2)
In [102]:
# Fit grid search
best_model = clf.fit(X_train, y_train)
In [103]:
print('Best hyper parameter:', best_model.best_params_)
In [104]:
# View best hyperparameters
training_score = best_model.best_score_
print(best_model.best_score_)
print(best_model.best_params_)
# Actual model object fit with those best parameters
# Shows default parameters that we did not specify
print(best_model.best_estimator_)
In [105]:
# Apply best params to fit this model and see training scores
#clfg = RandomForestRegressor(random_state=1, max_depth=10, n_estimators=100, max_features='auto',min_samples_leaf=1, min_samples_split=2)
clfg = best_model.best_estimator_
print("R Squared:",cross_val_score(clfg, X_train, y_train.values.ravel(), cv=10, scoring='r2').mean())
In [106]:
clfg = clfg.fit(X_train, y_train.values.ravel())
In [107]:
test_RMSE = np.sqrt(mean_squared_error(y_test, clfg.predict(X_test)))
print("Test RMSE: ", np.sqrt(mean_squared_error(y_test, clfg.predict(X_test))))
test_score = r2_score(y_test, clfg.predict(X_test))
print("Test R^2: ", r2_score(y_test, clfg.predict(X_test)))
In [108]:
#Store the accuracy results for each model in a dataframe for final comparison
r_df = pd.DataFrame({'Method':['Random Forest with RandomSearch CV'],'Training Score': training_score.round(3),'Testing score': test_score.round(3),
'RMSE': test_RMSE})
results_df = results_df.append(r_df)
results_df
Out[108]:
Hyperparameter Tuning - GridSearchCV - Gradient Boost¶
In [109]:
# Create the GridSearch estimator along with a parameter object containing the values to adjust
from sklearn.model_selection import GridSearchCV
n_estimators = [200, 300, 400, 500,1000]
learning_rate = [0.001,0.05,0.01,0.1]
max_depth = [8,10,20,30,40]
subsample = [0.5, 0.6, 0.8, 0.9, 1.0]
param_grid = dict(n_estimators = n_estimators, max_depth = max_depth,learning_rate = learning_rate,
subsample = subsample)
In [110]:
# Create grid search using 5-fold cross validation
clf = GridSearchCV(gboost, param_grid, cv=5, n_jobs = -1,verbose=2)
In [111]:
# Fit grid search
gboost_best_model = clf.fit(X_train, y_train)
In [112]:
# View best hyperparameters
training_score = gboost_best_model.best_score_
print(gboost_best_model.best_score_)
print(gboost_best_model.best_params_)
# Actual model object fit with those best parameters
# Shows default parameters that we did not specify
print(gboost_best_model.best_estimator_)
In [113]:
# Apply best params to fit this model and see training scores
#clfg = GradientBoostingRegressor(random_state=1, max_depth=5, n_estimators=500, max_features='sqrt',learning_rate=0.05,subsample=0.5)
clfg = gboost_best_model.best_estimator_
print("R Squared:",cross_val_score(clfg, X_train, y_train.values.ravel(), cv=10, scoring='r2').mean())
In [114]:
clfg = clfg.fit(X_train, y_train.values.ravel())
In [115]:
test_RMSE = np.sqrt(mean_squared_error(y_test, clfg.predict(X_test)))
print("Test RMSE: ", np.sqrt(mean_squared_error(y_test, clfg.predict(X_test))))
test_score = r2_score(y_test, clfg.predict(X_test))
print("Test R^2: ", r2_score(y_test, clfg.predict(X_test)))
In [116]:
#Store the accuracy results for each model in a dataframe for final comparison
r_df = pd.DataFrame({'Method':['Gradient Boost with GridSearch CV'],'Training Score': training_score.round(3),'Testing score': test_score.round(3),
'RMSE': test_RMSE})
results_df = results_df.append(r_df)
results_df
Out[116]:
Hyperparameter Tuning - RandomizedSearchCV - Gradient Boost¶
In [117]:
from sklearn.model_selection import RandomizedSearchCV
n_estimators = [200, 300, 400, 500,1000]
learning_rate = [0.001,0.05,0.01,0.1]
max_depth = [8,10,20,30,40]
subsample = [0.5, 0.6, 0.8, 0.9, 1.0]
param_grid = dict(n_estimators = n_estimators, max_depth = max_depth,learning_rate = learning_rate,
subsample = subsample)
In [118]:
# Create grid search using 5-fold cross validation
clf = RandomizedSearchCV(gboost, param_grid, cv=5, n_jobs = -1,verbose=2)
In [119]:
# Fit grid search
gboost_best_model = clf.fit(X_train, y_train)
In [120]:
# View best hyperparameters
training_score = gboost_best_model.best_score_
print(gboost_best_model.best_score_)
print(gboost_best_model.best_params_)
# Actual model object fit with those best parameters
# Shows default parameters that we did not specify
print(gboost_best_model.best_estimator_)
In [121]:
# Apply best params to fit this model and see training scores
#clfg = GradientBoostingRegressor(random_state=1, max_depth=6, n_estimators=200, max_features='sqrt',learning_rate=0.05,subsample=0.5)
clfg = gboost_best_model.best_estimator_
print("R Squared:",cross_val_score(clfg, X_train, y_train.values.ravel(), cv=10, scoring='r2').mean())
In [122]:
clfg = clfg.fit(X_train, y_train.values.ravel())
In [123]:
test_RMSE = np.sqrt(mean_squared_error(y_test, clfg.predict(X_test)))
print("Test RMSE: ", np.sqrt(mean_squared_error(y_test, clfg.predict(X_test))))
test_score = r2_score(y_test, clfg.predict(X_test))
print("Test R^2: ", r2_score(y_test, clfg.predict(X_test)))
In [124]:
#Store the accuracy results for each model in a dataframe for final comparison
r_df = pd.DataFrame({'Method':['Gradient Boost with RandomSearch CV'],'Training Score': training_score.round(3),'Testing score': test_score.round(3),
'RMSE': test_RMSE})
results_df = results_df.append(r_df)
results_df
Out[124]:
In [125]:
import xgboost as xgb
from xgboost.sklearn import XGBRegressor
#n_estimators = [50, 100, 150, 200]
#learning_rate = [0.05, 0.10, 0.15, 0.20, 0.25, 0.30 ]
#max_depth = [ 3, 4, 5, 6, 8, 10, 12, 15]
#min_child_weight = [ 1, 3, 5, 7 ]
#gamma = [ 0.0, 0.1, 0.2 , 0.3, 0.4 ]
#colsample_bytree = [ 0.3, 0.4, 0.5 , 0.7 ]
n_estimators = [50, 100, 150, 200]
learning_rate = [0.01, 0.1, 0.2, 0.3]
max_depth = range(3, 10)
colsample_bytree = [i/10.0 for i in range(1, 3)]
gamma = [i/10.0 for i in range(3)]
param_grid = dict(n_estimators = n_estimators, learning_rate = learning_rate, max_depth = max_depth,
colsample_bytree = colsample_bytree,
gamma = gamma)
In [126]:
# Create grid search using 5-fold cross validation
clf = GridSearchCV(xgboost, param_grid, cv=5, n_jobs = -1,verbose=2)
In [127]:
# Fit grid search
xgboost_best_model = clf.fit(X_train, y_train)
In [128]:
# View best hyperparameters
training_score = xgboost_best_model.best_score_
print(xgboost_best_model.best_score_)
print(xgboost_best_model.best_params_)
# Actual model object fit with those best parameters
# Shows default parameters that we did not specify
print(xgboost_best_model.best_estimator_)
In [129]:
# Apply best params to fit this model and see training scores
clfg = xgboost_best_model.best_estimator_
print("R Squared:",cross_val_score(clfg, X_train, y_train.values.ravel(), cv=10, scoring='r2').mean())
In [130]:
clfg = clfg.fit(X_train, y_train.values.ravel())
In [131]:
test_RMSE = np.sqrt(mean_squared_error(y_test, clfg.predict(X_test)))
print("Test RMSE: ", np.sqrt(mean_squared_error(y_test, clfg.predict(X_test))))
test_score = r2_score(y_test, clfg.predict(X_test))
print("Test R^2: ", r2_score(y_test, clfg.predict(X_test)))
In [132]:
#Store the accuracy results for each model in a dataframe for final comparison
r_df = pd.DataFrame({'Method':['XGBoost with GridSearch CV'],'Training Score': training_score.round(3),'Testing score': test_score.round(3),
'RMSE': test_RMSE})
results_df = results_df.append(r_df)
results_df
Out[132]:
Bootstrap Sampling¶
In [133]:
#Bootstrap Sampling
from sklearn.utils import resample
values = concrete_df_scaled.values
n_iterations = 500 # Number of bootstrap samples to create
n_size = int(len(concrete_df_scaled)) # size of a bootstrap sample
# run bootstrap
stats = list() # empty list that will hold the scores for each bootstrap iteration
for i in range(n_iterations):
# prepare train and test sets
train = resample(values, n_samples=n_size) # Sampling with replacement
# picking rest of the data not considered in sample
test = np.array([x for x in values if x.tolist() not in train.tolist()])
# fit model
gboost_boot = gboost_best_model.best_estimator_
gboost_boot.fit(train[:,:-1], train[:,-1]) # fit against independent variables and corresponding target values
y_test = test[:,-1] # Take the target column for all rows in test set
# evaluate model
predictions = gboost_boot.predict(test[:, :-1]) # predict based on independent variables in the test data
score = gboost_boot.score(test[:, :-1] , y_test)
stats.append(score)
In [134]:
# Plot confidence Interval
plt.figure(figsize=(20,5))
# plt.hist(stats)
sns.distplot(stats, bins=11)
alpha = 0.95 # for 95% confidence interval
p = ((1.0-alpha)/2.0) * 100 # tail regions on right and left .25 on each side indicated by P value (border)
lower = max(0.0, np.percentile(stats, p))
p = (alpha+((1.0-alpha)/2.0)) * 100
upper = min(1.0, np.percentile(stats, p))
plt.axvline(x=lower, color='r', linestyle=':')
plt.axvline(x=upper, color='r', linestyle=':', label='Confidence Interval Threshold')
plt.legend(loc="best",prop={"size":14})
print("With %.1f confidence interval, Gradient Boosting's score varies between %.1f%% and %.1f%%" %
(alpha*100, lower*100, upper*100))
With 100 samples, 95.0 confidence interval, Gradient Boosting's score varies between 88.4% and 93.1%
Model Building with Feature Extraction using PCA¶
In [135]:
from sklearn.decomposition import PCA
In [136]:
concrete_df_pca = concrete_df_scaled.copy()
concrete_df_pca
Out[136]:
In [137]:
# Initialize and fit pca
# Create a covariance matrix and calculate eigen values
from sklearn.decomposition import PCA
pca = PCA().fit(concrete_df_pca.drop('strength', axis=1))
In [138]:
# cumulative sum of variance explained with [n] features
eigen_vals = np.round(pca.explained_variance_ratio_, decimals=3)*100
np.cumsum(eigen_vals)
Out[138]:
In [139]:
print(pca.explained_variance_)
In [140]:
print(pca.components_)
In [141]:
# calculate variance ratios
var = pca.explained_variance_ratio_;var
Out[141]:
In [142]:
# Visualize the contribution of variance of each attribute
plt.figure(figsize=(10,8))
plt.bar(list(range(8)),pca.explained_variance_ratio_,alpha=0.5, align='center');
plt.ylabel('Variation')
plt.xlabel('Number of components')
Out[142]:
In [143]:
# 6 attributes explaines 95% of variance overall
plt.figure(figsize=(10,8))
plt.step(list(range(8)),np.cumsum(pca.explained_variance_ratio_), where='mid');
plt.ylabel('Cumulative of variation')
plt.xlabel('Number of components')
Out[143]:
There are 6 features which explains more than 95% of variance cumulatively in the dataset
In [144]:
# Create a new matrix using the n components
X = concrete_df_pca.drop('strength', axis=1)
X_proj = PCA(n_components=6).fit_transform(X)
y = concrete_df_pca['strength']
X_proj.shape
Out[144]:
In [145]:
# Divide the original and the projected dataset into 70:30 ration for train and test
X_proj_train, X_proj_test, y_train, y_test = train_test_split(X_proj, y, test_size=0.3, random_state=1)
Linear Regression¶
In [146]:
linreg = LinearRegression()
In [147]:
linreg_scores = cross_val_score(linreg, X_proj_train, y_train, cv=KFold(10, random_state = 1,shuffle = True),scoring='neg_mean_squared_error')
linreg.fit(X_proj_train, y_train)
print(linreg.coef_)
print()
# Calculate the score of Linear Regression
print('Training score :', linreg.score(X_proj_train, y_train))
print('Testing score :', linreg.score(X_proj_test, y_test))
# Calculate average RMSE
print('RMSE:',np.sqrt(-(linreg_scores)).mean())
Decision Tree with K-fold cross validation¶
In [148]:
dtree = DecisionTreeRegressor(random_state=1)
dtree_scores = cross_val_score(dtree, X_proj_train, y_train, cv=KFold(10, random_state = 1,shuffle = True), scoring='neg_mean_squared_error')
dtree.fit(X_proj_train, y_train)
Out[148]:
In [149]:
#Making the prediction
y_predict = dtree.predict(X_proj_test)
In [150]:
# Calculate the score of Decision Tree Regressor
print('Training Score:',dtree.score(X_proj_train,y_train))
print('Testing Score:',dtree.score(X_proj_test,y_test))
# calculate the average RMSE
print('RMSE:',np.sqrt(-(dtree_scores)).mean())
Random Forest with K-fold cross validation¶
In [151]:
rforest = RandomForestRegressor(random_state=1)
rforest_scores = cross_val_score(rforest, X_proj_train, y_train, cv=KFold(10, random_state = 1,shuffle = True), scoring='neg_mean_squared_error')
rforest = rforest.fit(X_proj_train, y_train)
In [152]:
y_predict = rforest.predict(X_proj_test)
In [153]:
# Calculate the score of Random Forest Regressor
print('Training Score:',rforest.score(X_proj_train,y_train))
print('Testing Score:',rforest.score(X_proj_test,y_test))
# calculate the average RMSE
print('RMSE:',np.sqrt(-(rforest_scores)).mean())
Gradient Boost with K-fold cross validation¶
In [154]:
gboost = GradientBoostingRegressor(random_state=1)
gboost_scores = cross_val_score(gboost, X_proj_train, y_train, cv=KFold(10, random_state = 1,shuffle = True), scoring='neg_mean_squared_error')
gboost = gboost.fit(X_proj_train, y_train)
In [155]:
y_predict = gboost.predict(X_proj_test)
In [156]:
# Calculate the score of Gradient Boost Regressor
print('Training Score:',gboost.score(X_proj_train,y_train))
print('Testing Score:',gboost.score(X_proj_test,y_test))
# calculate the average RMSE
print('RMSE:',np.sqrt(-(gboost_scores)).mean())
It looks like all the regressors tried above with PCA performed very poorly in terms of training and test scores and with high RMSE scores on PCA projected dataset and as such are able to predict the compressive strength of the cement with an acceptable accuracy.
Model Building with Feature Generation using Polynomial Features¶
In [157]:
from sklearn.preprocessing import PolynomialFeatures
In [158]:
concrete_df_poly = concrete_df_new.copy()
concrete_df_poly
Out[158]:
In [159]:
# Instantiate an poly object with 2nd degree
poly = PolynomialFeatures(degree=2, interaction_only=True)
In [160]:
# Fit and transform the X or input features
concrete_df_poly = poly.fit_transform(concrete_df_poly.drop('strength', axis=1))
In [161]:
# Get the shape of the newly generated features
print("Number of rows :",concrete_df_poly.shape[0])
print("Number of columns :",concrete_df_poly.shape[1])
In [162]:
# Join the strength column to create a polynomial dataset
concrete_df_poly = pd.DataFrame(concrete_df_poly,columns=poly.get_feature_names())
concrete_df_poly['strength'] = concrete_df_new['strength']
concrete_df_poly.shape
Out[162]:
In [163]:
# Visualize the correlation matrix
plt.figure(figsize=(15,10))
sns.heatmap(concrete_df_poly.corr())
Out[163]:
In [164]:
x = concrete_df_poly.corr()['strength']
upper_limit = 0.6
lower_limit = -0.5
fltr = list(map(lambda x: True if(x >= upper_limit or x <= lower_limit) else False, x))
x[fltr].sort_values()
Out[164]:
This polynomial features dataset (with a degree 2) has very less correlationship coming out among all the features and with the strength column indicating that no newly generated features are a very good strong predictor for compressive strength. All have a coefficient less than or close to 0.5.
Model Building with Automatic Feature Engineering¶
In [165]:
# Import package
import featuretools as ft
In [166]:
ft.list_primitives().tail(30)
Out[166]:
In [167]:
# initialize entityset
es = ft.EntitySet(id = 'concrete_entitySet')
In [168]:
concrete_df_new
Out[168]:
In [171]:
# initialize entity
es.entity_from_dataframe(entity_id = 'concrete_entity_id', index = 'id',dataframe = concrete_df_new, make_index = True)
Out[171]:
In [172]:
es['concrete_entity_id'].df
Out[172]:
In [173]:
es['concrete_entity_id']
Out[173]:
In [174]:
concrete_new_features, feature_defs = ft.dfs(entityset = es, target_entity = 'concrete_entity_id',
trans_primitives = ['add_numeric_scalar', 'multiply_numeric','divide_numeric'])
concrete_new_features.head()
Out[174]:
In [175]:
# Visualize the correlation matrix
plt.figure(figsize=(15,10))
sns.heatmap(concrete_new_features.corr())
Out[175]:
In [176]:
# Correlation matrix
corr = concrete_new_features.corr()
corr
Out[176]:
In [177]:
concrete_new_features
Out[177]:
In [178]:
x = concrete_new_features.corr()['strength']
upper_limit = 0.6
lower_limit = -0.5
fltr = list(map(lambda x: True if(x >= upper_limit or x <= lower_limit) else False, x))
x[fltr].sort_values()
Out[178]:
Extracting the composite features from above which will work with strength and adding them to the dataset
In [179]:
cols = list(concrete_df_new.columns)
cols.remove('id')
cols
Out[179]:
In [180]:
cols.append('water / age')
cols.append('fineagg / age')
cols.append('coarseagg / age')
print(cols)
In [181]:
concrete_ft_eng = concrete_new_features[cols]
concrete_ft_eng
Out[181]:
In [182]:
std_scale = StandardScaler()
concrete_ft_scaled = pd.DataFrame(std_scale.fit_transform(concrete_ft_eng),columns = concrete_ft_eng.columns)
In [183]:
X = concrete_ft_scaled.drop(['strength'], axis=1)
Y = concrete_ft_scaled[['strength']]
In [184]:
# Split into training and test set
from sklearn.model_selection import train_test_split
X_train, X_test, y_train, y_test = train_test_split(X, Y, test_size=0.3,random_state=1)
Random Forest with K-fold cross validation¶
In [185]:
rforest = RandomForestRegressor(random_state=1)
rforest_scores = cross_val_score(rforest, X_train, y_train, cv=KFold(10, random_state = 1,shuffle = True), scoring='neg_mean_squared_error')
rforest = rforest.fit(X_train, y_train)
In [186]:
y_predict = rforest.predict(X_test)
# Calculate the score of Random Forest Regressor
print('Training Score:',rforest.score(X_train,y_train))
print('Testing Score:',rforest.score(X_test,y_test))
# calculate the average RMSE
print('RMSE:',np.sqrt(-(rforest_scores)).mean())
Gradient Boost with K-fold cross validation¶
In [187]:
gboost = GradientBoostingRegressor(random_state=1)
gboost_scores = cross_val_score(gboost, X_train, y_train, cv=KFold(10, random_state = 1,shuffle = True), scoring='neg_mean_squared_error')
gboost = gboost.fit(X_train, y_train)
In [188]:
y_predict = gboost.predict(X_test)
In [189]:
# Calculate the score of Gradient Boost Regressor
print('Training Score:',gboost.score(X_train,y_train))
print('Testing Score:',gboost.score(X_test,y_test))
# calculate the average RMSE
print('RMSE:',np.sqrt(-(gboost_scores)).mean())
Random forest model seems to be highly overfit as there is quite a drop in the accuracy of the test set and Gradientboosregressor performs better than randomforst regressor and hence it can be further tuned to squeeze the extra performance
Hyperparameter Tuning - GridSearchCV - Gradient Boost¶
In [190]:
n_estimators = [200, 300, 400, 500,1000]
learning_rate = [0.001,0.05,0.01,0.1]
max_depth = [8,10,20,30,40]
subsample = [0.5, 0.6, 0.8, 0.9, 1.0]
In [191]:
param_grid = dict(n_estimators = n_estimators, max_depth = max_depth,learning_rate = learning_rate,
subsample = subsample)
In [192]:
# Create grid search using 5-fold cross validation
clf = GridSearchCV(gboost, param_grid, cv=5, n_jobs = -1,verbose=2)
In [193]:
# Fit grid search
gboost_best_model = clf.fit(X_train, y_train)
In [194]:
# View best hyperparameters
training_score = gboost_best_model.best_score_
print(gboost_best_model.best_score_)
print(gboost_best_model.best_params_)
In [195]:
# Actual model object fit with those best parameters
# Shows default parameters that we did not specify
print(gboost_best_model.best_estimator_)
In [196]:
# Apply best params to fit this model and see training scores
clfg = gboost_best_model.best_estimator_
print("R Squared:",cross_val_score(clfg, X_train, y_train.values.ravel(), cv=10, scoring='r2').mean())
In [197]:
clfg = clfg.fit(X_train, y_train.values.ravel())
In [198]:
test_RMSE = np.sqrt(mean_squared_error(y_test, clfg.predict(X_test)))
print("Test RMSE: ", np.sqrt(mean_squared_error(y_test, clfg.predict(X_test))))
test_score = r2_score(y_test, clfg.predict(X_test))
print("Test R^2: ", r2_score(y_test, clfg.predict(X_test)))
In [199]:
#Bootstrap Sampling
from sklearn.utils import resample
values = concrete_ft_scaled.values
n_iterations = 500 # Number of bootstrap samples to create
n_size = int(len(concrete_ft_scaled)) # size of a bootstrap sample
# run bootstrap
stats = list() # empty list that will hold the scores for each bootstrap iteration
for i in range(n_iterations):
# prepare train and test sets
train = resample(values, n_samples=n_size) # Sampling with replacement
# picking rest of the data not considered in sample
test = np.array([x for x in values if x.tolist() not in train.tolist()])
# fit model
gboost_boot = gboost_best_model.best_estimator_
gboost_boot.fit(train[:,:-1], train[:,-1]) # fit against independent variables and corresponding target values
y_test = test[:,-1] # Take the target column for all rows in test set
# evaluate model
predictions = gboost_boot.predict(test[:, :-1]) # predict based on independent variables in the test data
score = gboost_boot.score(test[:, :-1] , y_test)
stats.append(score)
In [200]:
# Plot confidence Interval
plt.figure(figsize=(20,5))
# plt.hist(stats)
sns.distplot(stats, bins=11)
alpha = 0.95 # for 95% confidence interval
p = ((1.0-alpha)/2.0) * 100 # tail regions on right and left .25 on each side indicated by P value (border)
lower = max(0.0, np.percentile(stats, p))
p = (alpha+((1.0-alpha)/2.0)) * 100
upper = min(1.0, np.percentile(stats, p))
plt.axvline(x=lower, color='r', linestyle=':')
plt.axvline(x=upper, color='r', linestyle=':', label='Confidence Interval Threshold')
plt.legend(loc="best",prop={"size":14})
print("With %.1f confidence interval, Gradient Boosting's score varies between %.1f%% and %.1f%%" %
(alpha*100, lower*100, upper*100))
In [ ]: